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Abstract 

Using molecular dynamic simulations we study three families of continuous core-softened poten- 
tials consisting of two length scales: a shoulder scale and an attractive scale. All the families have 
the same slope between the two length scales but exhibit different potential energy gap between 
them. For each family three shoulder depths are analyzed. We show that all these systems exhibit 
a liquid-liquid phase transition between a high density liquid phase and a low density liquid phase 
ending at a critical point. The critical temperature is the same for all cases suggesting that the 
critical temperature is only dependent on the slope between the two scales. The critical pressure 
decreases with the decrease of the potential energy gap between the two scales suggesting that the 
pressure is responsible for forming the high density liquid. We also show, using the radial distri- 
bution function and the excess entropy analysis, that the density, the diffusion and the structural 
anomalies are present if particles move from the attractive scale to the shoulder scale with the 
increase of the temperature indicating that the anomalous behavior depends only in what happens 
up to the second coordination shell. 

PACS numbers: 64.70.Pf, 82.70.Dd, 83.10.Rs, 61.20.Ja 
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I. INTRODUCTION 



The phase behavior of single component systems as particles interacting via the so-called 
core-softened (CS) potentials is receiving a lot of attention recently. These potentials exhibit 
a repulsive core with a softening region with a shoulder or a rampi - -. These models originates 
from the desire of constructing a simple two-body isotropic potential capable of describing 
the complicated features of systems interacting via anisotropic potentials. This procedure 
generates systems that are analytically and computationally tractable and that one hopes 
are capable to retain the qualitative features of the real complex systems^— . 

The physical motivation behind these studies is the recently acknowledged possibility 
that some single component systems interacting through a core-softened potential display 
density and diffusion anomalies. This opened the discussion about the relation between the 
existence of thermodynamic anomalies in liquids and the form of the effective potential^ - —. 

These anomalies appear in two different ways. First, it is the density anomaly. Most 
liquids contract upon cooling. This is not the case of water and other fluid systems. For 
water the specific volume at ambient pressure starts to increase when cooled below T « 4°C. 
The anomalous behavior of water was first suggested 300 years ago^ and was confirmed by 
a number of experiments^ 1 ^. Besides, between 0.1 MPa and 190 MPa water also exhibits 
an anomalous increase of compressibility^ 1 ^ and, at atmospheric pressure, an increase of 
isobaric heat capacity upon cooling^ 1 ^ 3 -. 

Experiments for Te™^ Ga, Bi,— S,— »2Z and Ge^Tess,— and simulations for silica,— 1 ^— 
silicon 31 and BeF 2 ,— show the same density anomaly. 

But density anomaly is not the only unusual behavior that these materials have. For a 
normal liquid the diffusion constant, D, decreases under compression. This is not the case 
of water. D increases on compression at low temperature, T, up to a maximum D max (T) 
at p — PDma.x{T)-^22.- Numerical simulations for SPC/E water— recover the experimental 
results and show that the anomalous behavior of D extends to the metastable liquid phase 
of water at negative pressure, a region that is difficult to access for experiments^— . In this 
region the diffusivity D decreases for decreasing p until it reaches a minimum value D min (T) 
at some pressure pDmm{T), and the normal behavior, with D increasing for decreasing p, is 
reestablished only for p < pDmm{T)^~—. Besides water, silica^ 1 ^ and silicon^ also exhibit 
a diffusion anomalous region. 
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Acknowledging that CS potentials may engender density and diffusion anomalous behav- 
ior, a number of CS potentials were proposed to model the anisotropic systems described 
above. They possess a repulsive core that exhibits a region of softening where the slope 
changes dramatically This region can be a shoulder or a ramp-&i£»i&^~— . These models 
exhibit density and diffusion anomalies, but depending on the specific shape of the poten- 
tial, the anomalies might be hidden in the metastable phased. Also there are a number of 
core-softened potentials in which the anomalies are not present§2^. The relation between 
the specific shape of the effective core-softened potential and the microscopic mechanism 
necessary for the presence of the anomalies is still under debate^ 1 ^ 1 ^ 1 ^. 

Recently it was suggested that the link between the presence of the density and diffusion 
anomaly and the microscopic details of the system can be analyzed in the framework of 
the excess-entropy-based formalism 64 applied to similar systems by Errington et al— and 
Chakraborty and Chakravarty^. Within this approach the presence of the density and the 
diffusion anomalies are related to the density dependence of the excess entropy, s ex . 

The computation of the excess entropy, however, requires integrating the radial distri- 
bution function in the whole space. The anomalous behavior, however, seems to depend 
only in the two length scales present in the system and, therefore, should not depend on 
the particle distributions far away. Here we propose that two length scales potentials will 
have density and diffusion anomalies if the two scales would be accessible. In principle the 
accessibility only depends in the distribution of particles in these two distances. Therefore, 
the knowledge of the complete excess entropy is not necessary for knowing if the system 
has or not anomalies. The behavior of the partial excess entropy, computed only up to 
the second coordination shell should give enough information to determine if a system has 
anomalies or not. 

In this paper we test this assumption by computing the pressure-temperature phase 
diagram and the excess entropy for three families of core-softened potentials that have two 
length scale: a shoulder scale and an attractive scale^ 1 ^. In all the three families the 
slope between the two scales is the sameJ^. The shoulder scale is made more favorable 
by decreasing the energy gap between the two length scales (see potentials A, B and C in 
Fig. [1]). In addition the shoulder scale becomes more favorable by making the depth of the 
shoulder scale deeper (see potentials Al, A2 and A3 in Fig. [p. The slope between the two 
length scales is kept fixed in order to have the liquid-liquid critical point and the density 
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anomalous region in the same region of the pressure temperature phase diagram^. 

The remaining of this paper goes as follows. In Sec. [Ill the model is introduced. The 
simulations details are given Sec. IIHI In Sec. [IV] the results are discussed. Finally, Sec. [V] 
presents the conclusion. 



II. THE MODEL 



We consider a system of iV particles, with diameter a, where the pair interaction is 
described by a family of continuous potentials given by 
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(a)Values of ci = 0.90c re/ , 
hi = 0.30/T e/ ,0.60/T e/ and 
0.90/i re/ . 



(c) Values of ci = 1.10c re ^, 
hi = 0.22/r e /, 0-44h re/ and 
0Mh re f. 



(b) Values of ci = 1.00c re/ , 
hi = 0.25/i re/ ,0.50/i re/ and 
0.75h re f. 

FIG. 1: Interaction potential obtained by changing parameters hi in the Eq. (TjJ.The potential 
and the distances are in dimensionless units U* = C//7 and r* = r/ro-Here we use e/7 = 0.02 and 
cr/r = 1.47. 

The first term is a Lennard- Jones potential and the second term is composed by four 
Gaussians, each Gaussian is centered in Cj. This potential can represent a whole family 
of intermolecular interactions, depending of the choice of the parameters a, b, a, {hj, Cj,Wj}, 
with j — 1, ... ,4. The values of these parameters are in Table [B The parameters are chosen 
in order to obtain a two length scale potential^ - — that is related to the interaction between 
two tetramers^ 1 ^. 
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(a) Distance between the scales. (b) Different depths of the first scale. 

FIG. 2: Scheme to distinguish the distance between the scales and depths of the first scale. 

The simulations are made in dimensionless units, therefore all the physical quantities are 
given in terms of the energy scale 7 and the distance scale where 7 is the energy scale 
and r is the length scale chosen so the minimum of the potential in the B 3 case is about 
r* = 1. The distance scale is chosen so the minimum of the shoulder scale in the B 3 case 
Here we use e/7 = 0.02 and cr/r = 1.47. 

Modifying the parameters c\ and h\ in the Eq. ([[]), according to Table ILT1 allow us to 
change the depth of the shoulder well, as illustrated in Fig. [TJ Here we use nine different 
values for hi and they are expressed as a multiple of a reference value h r ^ . We also use 
three different values of C\ and they are expressed as a multiple of a reference value . 
For all the nine cases the values of a, b, {wj} with j — 1, . . . , 4, q with i — 2, ... ,3, h re f and 
c hre f . Table U gives the parameter values in A and kcal/mol consistent with modeling ST4 
water—. 

Modifying the distance between the two minima of the two scales, shoulder scale (SC) 
and attractive scale (AS), leads to the three families A, B and C as shown in Fig. [TJ The 
changes in the distance between the two length scales were done in such way to preserve the 
slope between the two scales and, therefore, to have in all the cases the region of density 
anomaly in the same region of the pressure-temperature phase diagram as proposed by Yan 
et. aP. 

The family A has the largest distance and the largest potential energy gap between the 
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two length scales (a a = 0.72 in Fig. [2]), the family B has the intermediate distance and 
intermediate potential energy gap between the two length scales (ob = 0.62 in Fig. H}, and 
the family C has the shortest distance and the smallest potential energy gap between the 
two length scales, (ac = 0.52 in Fig. [2]). 

For each family, we analyze three different depths of the shoulder scale represented by 1 , 
2 and 3 (see Fig. [2]). The potentials 1 have the most shallow shoulder scale , the potentials 
2 have intermediate shoulder depth and the potentials 3 have the deepest shoulder scale. 
Table ITO1 gives the values for the depths for each one of the families. 

In summary, we analyze nine different potentials: A±, A 2 , A 3 , Bi, B 2 , B 3 , Cj, C 2 and C 3 as 
illustrated in Fig. [TJ The values of the different hi for each case are listed in Table [TTJ 

Barraz et al. 18 investigated the family B. It was shown that this potential exhibits 
thermodynamic, dynamic and structural anomalies if the shoulder scale is not too deep. 
Their result suggests that in order to have anomalies it is necessary but no sufficient to have 
two length scales competing. Both length scales must be accessible. When the shoulder 
scale becomes too deep the particles are trapped in this length scale and no anomaly is 
present. By making the shoulder deeper we are decreasing the difference in energy between 
the scales and therefore destroying the competition. 

Here we explore the accessibility or the absence of accessibility by changing both the 
energy gap between the two length scales and the distance between them. 

III. THE SIMULATION DETAILS 

The properties of the system were obtained by NVT molecular dynamics using Nose- 
Hoover heat-bath with coupling parameter Q — 2. The system is characterized by 500 
particles in a cubic box with periodic boundary conditions, interacting with the intermolec- 
ular potential described above. All physical quantities are expressed in reduced units^ given 
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by 



t* = tr (m/ 7 ) 1/2 
rp* _ k B T 

7 

p * = PA 

7 

* 3 

P = pr 



D* = D x — 2 . (2) 
V 7^o 

Standard periodic boundary conditions together with predictor-corrector algorithm were 
used to integrate the equations of motion with a time step At* = 0.002 and potential cut 
off radius r* = 3.5. The initial configuration were set on solid and on liquid states and, 
in both cases, the equilibrium state was reached after t* eq = 1000 (what is in fact 500, 000 
steps since At* = 0.002). From this time on the physical quantities were stored in intervals 
of At* R = 1 during t* R = 1000. The system is uncorrelated after t* d = 10, from the velocity 
auto-correlation function. With 50 descorrelated samples were used to get the average of 
the physical quantities. 

The thermodynamic stability of the system was checked by analyzing the dependence of 
pressure on density namely 

|<° < 3 > 

, by the behavior of the energy and also by visual analysis of the final structure, searching 
for cavitation. 

At the phase boundary between the liquid and the amorphous phases we found stable 
states points at both phases. The state with the lower energy was considered. In this partic- 
ular region of the pressure-temperature phase diagram the energy was a good approximation 
for the Helmoltz free energy. 

The error bars for temperatures and pressures away from the critical region are smaller 
than the size of the gray lines. The error bar near the critical point are AT = 0.0025 and 
Ap = 0.05. Our error is controlled by making averages of uncorrelated measures. 
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TABLE I: Parameters for potentials A, B and C in reduced unit. 



Parameter Value 


Parameter Value 


Parameter Value 


Parameter Value 


a 9.056 
b 4.044 
e 0.020 
a 1.475 


wi 0.085 
w 2 0.618 
w 3 0.826 
Wi 0.214 


c[ ef 0.996 
c 2 0.529 
c 3 1.598 
c 4 1.929 


-3.79 
h 2 1.209 
/i 3 -1.503 
h 4 0.767 



TABLE II: Parameters of c\ and /ii for potentials A, B and C. 



Potentials 


Values 


Potentials 


Values 


Potentials 


Values 


Potentials 


Values 


A 


0.90 c[ ef 


A x 


0.30 h[ ef 


Bi 


0.25 h[ ef 


Ci 


0.22 h[ ef 


B 


1.00 c[ ef 


A 2 


0.60 h[ ef 


B 2 


0.50 h\ ef 


c 2 


0.44 h\ e f 


C 


1.10 c[ ef 


A 3 


0.90 h[ ef 


B 3 


0.75 h\ e} 


C 3 


0.66 /ii e/ 



IV. RESULTS 

Pressure-Temperature Phase Diagram 

First, we explore the effects that the increase of the shoulder depth and the decrease of 
the distance between the two scales have in the location in the pressure-temperature phase 
diagram of the different phases. Fig. |3] illustrates the pressure versus temperature phase 
diagram of the three families A, B and C of potentials. At high temperatures there are a 
fluid phase and a gas phase (not shown). These two phases coexist at a first order line that 
ends at a critical point (see Table HVl for the pressure and the temperature values). 

At low temperatures and high pressures there are two liquid phases coexisting at a first 
order line ending at a second critical point (see Table IVl for the pressure and the temperature 



TABLE III: Values for the depths for each one of the families. 



Potential 


Value 


Potential 


Value 


Potential 


Value 




7.10 




4.94 


b Cl 


2.95 


bA 2 


6.20 


bB 2 


4.07 


bc 2 


2.33 


bA 3 


5.28 


bB 3 


3.32 


bc 3 


1.60 



s 



TABLE IV: First critical point location for potentials A, B and C. 



Potential T*x p* cl 


Potential T* ± p* cl 


Potential T* ± p* cl 


A x 1.94 0.074 
A 2 1.95 0.074 
A 3 1.97 0.076 


Bx 1.93 0.072 
B 2 1.98 0.078 
B 3 2.02 0.080 


Cx 1.98 0.076 
C 2 2.08 0.088 
C 3 2.20 0.099 


TABLE V: Second critical point location for potentials A, B and C. 


Potential T* 2 p* c2 


Potential T* 2 p* 2 


Potential T* 2 p* c2 


Ax 0.34 4.21 
A 2 0.48 3.16 
A 3 0.59 1.70 


Bx 0.35 3.44 
B 2 0.48 1.86 
B 3 0.57 0.49 


Cx 0.34 2.59 
C 2 0.41 0.78 
C 3 0.53 -0.63 



values). The thermodynamic stability of the state points was checked by analyzing the 
dependence of pressure on density using Eq. [3j The critical point was identified in the graph 
by the region where isochores cross. The coexistence line was obtained as the medium line 
between the stability limit of each phase. 

Comparison between the cases Ai, Bx and Cx indicates that as the distance between the 
two scales decreases, the critical pressure decreases but the critical temperature remains the 
same as illustrated in Fig. HJ This observation is confirmed in the cases A 2 , B 2 and C 2 and 
in the cases A 3 , B 3 and C 3 (see Table IVl for the critical pressures and the temperatures). 

How can this result be understood? The two liquid phases are formed due to the presence 
of the two competing scales. The low density liquid is related to the attractive scale while 
the high density liquid is related to the shoulder scale. In order to reach the high density 
liquid phase the system has to overcome a large potential energy but also have to become 
very compact in the case A\. The potential energy gap and the distance between the two 
length scales are higher in the case A than in the case C, therefore the pressure needed for 
forming a high density liquid phase is higher in the case A than it is in the case C. 
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FIG. 3: Pressure-temperature phase diagram for the family A in the left-hand side, for the family 
B in the middle and for C in the right-hand side. The thin solid lines are the isochores 0.30 < p* < 
0.65. The liquid-liquid critical point is the dot, the temperature of maximum density is the solid 
thick line, the diffusion extrema is the dashed line and the structural extrema is the dashed-dotted 
line. The dotted line indicates the limit between the fluid and the amorphous regions. For the 
diffusion coefficient and radial distribution functions it is possible to determine the amorphous 
region. The full line on the left side of the second critical point is an approximate line of the 
coexistence of fluids of high and low density. 
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TABLE VI: Values of pressure location in the amorphous region for cases A, B and C. 



Potentials 


Values 


Potentials 


Values 


Potentials 


Values 


Ax 


-0.90 < p* < 3.70 


Bi 


-0-91 <P*< 3.40 


Ci 


-0.85 <p*< 2.35 


A 2 


-0.30 < p* < 2.80 


B 2 


-0-89 <p*< 1.80 


c 2 


-0.76 <p*< 0.55 


A 3 


-1-05 <p*< 1.53 


B 3 


-1-00 <p*< 0.48 


c 3 
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FIG. 4: Location of the critical points on pressure-temperature phase diagram for cases A, B and 
C. 



At very low temperatures the system becomes less diffusive and crystallization might 
be expected. Here we do not explore all the possible crystal phases, but instead as the 
temperature is decreased from the liquid phase an amorphous phase is formed. The dotted 
line in Fig. |3] shows the separation between the fluid phase and the amorphous region. The 
amorphous phase is identified by the region where the diffusion coefficient becomes zero and 
the radial distribution functions do not exhibit the periodicity of a solid. Table [VI] have 
the characteristic pressure and temperature values of the amorphous phase boundary for 
the different models. It shows that the region in the pressure-temperature phase diagram 
where the amorphous phase is present shrinks as the shoulder part of the potential becomes 
deeper. 



Density Anomaly 

Next, we test the effects that the decrease of the distance between the two length scales 
and the increase of the depth of the first scale have in the presence of the density anomaly. 
Fig. [3] shows the isochores 0.30 < p* < 0.65 represented by thin solid lines for the nine 
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models. Equation 

indicates that the temperature of minimum pressure at constant density is the temperature 
of maximum density at constant pressure, the TMD. The TMD is the boundary of the region 
of thermodynamic anomaly, where a decrease in the temperature at constant pressure implies 
an anomalous increase in the density and therefore an anomalous behavior of density (similar 
to what happens in water). Figs. Elshow the TMD as solid thick lines. For all potentials A, 
B and the potentials C\ and C2 the TMD is present but for potential C3 the TMD is not 
observed. 

Similarly to what happens with the location of amorphous region, as the distance between 
the shoulder and the attractive scales decreases, the region in the pressure-temperature phase 
diagram delimited by the TMD line shrinks and disappears for the case C3. For the potential 
B the TMD line is located at temperatures bellow the temperature of the liquid-liquid critical 
point. The thermodynamic parameters that limits the TMD in phase diagram are shown 
in the Table IVTa where pi represents the values of (p*,T*,p*) for the point of the lowest 
pressure in the TMD line, p m is the point with the highest temperature and ph is the point 
with the highest pressure. 

How can this result be understood? The density anomalous behavior arises from the 
competition between the two length scales: the shoulder scale and the attractive scale. At 
high pressures the shoulder scale wins and at low pressures the attractive scale wins. The 
density anomalous region exists only in the intermediate pressure range where clusters of 
both scales are present. The value of the "high" pressure and of the "low" pressure is 
determined by the difference in energy between the two scales. If the difference is too small 
the low and high pressures are too close and no anomaly appears. 

Diffusion anomaly 

Then, we check the effects that the decrease of the distance between the two scales have 
in the location in the pressure-temperature phase diagram of the diffusion anomaly. The 
diffusion coefficient is obtained from the expression: 

D = lim (fc(* + t ) -^ifo)] 2 )*. ( 5 ) 
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TABLE VII: Limiting values for density (p*), temperature (T*) and pressure (p*) of the thermo- 
dynamics anomalies on pressure-temperature diagram. Here the point p\ represents the density, 
temperature and pressure of the point of the lowest pressure in the TMD line, p m represents the 
point of the highest temperature and ph represents the point of the highest pressure of the TMD 
line. 



cases pi p m p h 


cases pi p m p h 


cases pi p m p h 


p* 0.47 0.54 0.61 
A 1 T* 0.71 1.00 0.72 
p* 1.74 3.22 4.50 


p* 0.47 0.52 0.57 
B x T* 0.71 0.85 0.69 
p* 1.50 2.50 3.30 


p* 0.45 0.51 0.55 
Ci T* 0.56 0.69 0.40 
p* 0.79 1.96 2.60 


p* 0.46 0.51 0.59 
A 2 T* 0.70 1.00 0.54 
p* 1.36 2.20 3.17 


p* 0.46 0.50 0.54 
B 2 T* 0.67 0.76 0.63 
p* 0.90 1.40 1.80 


p* 0.46 0.48 0.50 
C 2 T* 0.49 0.53 0.48 
p* 0.23 0.45 0.64 


p* 0.45 0.48 0.54 
A 3 T* 0.66 0.89 0.69 
p* 0.91 1.27 1.72 


p* 0.40 0.42 0.43 
B 3 T* 0.44 0.54 0.52 
p* 0.15 0.29 0.36 


p* - - - 

C 3 T* - 

* 



where fj(t) are the coordinates of particle j at time t, and (• • • )t denotes an average over 
all particles and over all to- 

Fig. |5] shows the behavior of the dimensionless translational diffusion coefficient, D*, as 
function of the dimensionless density, p*, at constant temperature for the cases. A 3 , B 3 and 
C 3 . The solid lines are a polynomial fits to the data obtained by simulation (the dots in 
the Fig. [5]). For normal liquids, the diffusion at constant temperature increases with the 
decrease of the density. For the cases A\, A 2 , B\, B 2 , B 3 , C\, C 2 (not shown) and for the cases 
A 3 and B 3 the diffusion has a region in the pressure-temperature phase diagram in which 
the diffusion increases with density. This is the diffusion anomalous region illustrated in 
Fig. [3] as a dashed line. 

Similarly to what happens with the location of the TMD, as the two length scales becomes 
closer, the region in the pressure-temperature phase diagram delimited by the extrema of 
the diffusion goes to lower pressures, shrinks and disappears for the case C 3 . 

Fig. |3] illustrates that the region in the pressure-temperature phase diagram where the 
dynamic anomaly occurs englobes the region where the thermodynamic anomaly is present. 
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FIG. 5: Diffusion coefficient versus reduced density for T* = 0.6,0.7,0.8,0.9 from the bottom to 
the top for the case A3, T* = 0.4,0.5,0.6,0.7 from the bottom to the top for the case -B3 and 
T* = 0.5, 0.6, 0.7, 0.8, 0.9 from the bottom to the top for the case C3. The dots are the simulational 
data and the solid lines are polynomial fits. The dashed lines connect the densities of minima and 
maxima diffusivity that limit the diffusion anomalous region. 

This hierarchy between the anomalies is observed in simulations M£H£1 and in experiments^ 
for bulk water. 



Structural anomaly 

Now, we test the effects that the decrease of the distance between the two length scales 
have in the location in the pressure-temperature phase diagram of the structural anomalous 
region. 

The translational order parameter is denned g^^kJZ 

t= [ C \g(ti)-l\dti (6) 

JO 

where £ = rp5 is the distance r in units of the mean interparticle separation p~^, £ c is the 
cutoff distance set to half of the simulation box times^ p~^, g(£) is the radial distribution 
function which is proportional to the probability of finding a particle at a distance £ from a 
referent particle. The translational order parameter measures how structured is the system. 
For an ideal gas it is g — 1 and t — 0, while for the crystal phase it is g ^ 1 over long 
distances resulting in a large t. Therefore for normal fluids t increases with the increase of 
the density. 
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FIG. 6: The translational order parameter as a function of density for fixed temperatures: T* = 
1.50, 1.40, 1.30, 1.20, 1.10, 1.00, 0.90, 0.80, 0.70 and 0.60 (from top to bottom) for the setting 
potentials A and C; and T* = 1.00, 0.90, 0.80, 0.70 and 0.60 for the setting potential B. The 
dot-dashed lines locate the density of maxima e minima t*. 

Fig. [6] shows the translational order parameter as a function of the density for fixed 
temperatures for potentials A 3 , B 3 and C3. The dots represent the simulation data and 
the solid line the polynomial fit to the data. For the potentials A3 and B% there are a 
region of densities in which the translational parameter decreases as the density increases. 
A dotted-dashed line illustrates the region of local maximum and minimum of t* limiting 
the anomalous region. For the potential C3, t* increases with the density. No anomalous 
behavior is observed. The potentials A\, A 2 , Bi, B 2 , C\ and C 2 that do show anomalous 
behavior are not shown here for simplicity. 

Fig. [3] shows the structural anomaly for the families A, B and C, as dotted-dashed 
lines. It is observed that the region of structural anomaly embraces both dynamic and 
thermodynamic anomalies. As the distance between the shoulder and the attractive scales 
is decreased, the structural anomalous region in the pressure-temperature phase diagram 
shrinks. 

Another measure of the anomalous behavior is the orientational order parameter—, Qq. 
This parameter is used to get information about tetrahedral order of the molecules. For two 
length scales spherical symmetric continuous (continuous force) potentials, Q^-^ exhibit 
a region of temperatures in which it decreases with increasing density. The maximum of Q$ 
is located in the same region in the pressure temperature phase diagram of the maximum 
of the translational order parameter 55 . A similar behavior is expected for our potential. 
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In resume, for all the density, diffusion and structural anomalous regions in the pressure- 
temperature phase diagram, as the two length scales becomes closer the region of the 
pressure-temperature phase diagram occupied by the anomalous region shrinks. The same 
effect is observed when the shoulder scale becomes deeper—. 

Radial distribution function 

What is the origin of the disappearance of the thermodynamic, dynamic and structural 
anomalous behavior with the decrease of the distance between the scales? In order to 
answer to this question the behavior of the radial distribution function for the nine different 
potentials is studied. The radial distribution function is a measure of the probability of 
finding a pair of atoms separated by r. This function is defined as 



where fj(t) are the coordinates of particle % and j at time t, V is volume of system, N is 
number of particles and (•■■}* denotes an average over all particles. 

Recently it was shown that a necessary condition for the presence of density anomaly is 
to have particles moving from one scale to the other as the temperature is increased^^ 1 ^, 
for a fixed density. Here we test if this assumption is confirmed in the potentials we are 
analyzing. Fig. [7] illustrates the radial distribution function versus distance for a fixed 
density and various temperatures for the potentials A 3 , B 3 and C3. For the potentials A 3 
and B 3 the percentage of particles in the first length scale increases while the percentage of 
particles in the second length scale decreases as the temperature is increased. This means 
that as the system is heated at constant density particles move from one scale to the other. 
This behavior is also observed for the potentials A±, A2, Bi, B2, C% and C2 (not shown here 
for simplicity) and confirms our assumption that the presence of anomalies is related with 
particles moving from one length scale to the other length scale^ 1 ^. This is not the case for 
the potential C3. For the case C 3 , as the temperature is increased the particles move from 
the second to the other further away coordination shells and the percentage of particles at 
the first scale is not affected by the increase in temperature and therefore no anomaly is 
observed^ 1 ^. 
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How can we understand this result? The density anomaly appears if particles move from 
the second length scale to the first length scale. In the case of the potentials A3 and B3 the 
difference in energy between the two scales is high and heat is required for having particles 
reaching the shoulder length scale. Consequently, as the temperature is increased at constant 
density, more particles will be at the first scale and pressure decreases (see Fig. [7]). In the 
potential C3 the difference in energy between the two length scales is small. Almost no heat 
is required to have particles in the first length scales that saturates. So particles actually do 
not move from one scale to the other. 

This picture in terms of the presence of particles in the different shells can also be checked 
in the framework of the excess entropy^ 1 ^. 




FIG. 7: Radial distribution function versus reduced distance for three cases, A3, B3, C3 for reduced 
density p* = 0.480,0.420 and 0.390, respectively. In cases A3 and B3 the first peak of g(r*) 
increases with increasing temperature, while the second peak decreases. In C3 case for the first 
peak potential keeps constant while the second peak decreases. 



Excess entropy and anomalies 

The excess entropy is defined as the difference between the entropy of the real fluid and 
that of an ideal gas at the same temperature and density. It also can be given by its two 
body contribution of s ex , 

POO 

Ssex ~ s 2 = -2vrp / [g(r) In g(r) - g(r) + l}r 2 dr, (8) 
Jo 

gives a good approximation of s sex . 
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What can we learn from the excess entropy about the mechanism responsible for the 
density, the diffusion and the structural anomalies? In the sec. |IV]we have shown using the 
radial distribution function that for potentials that exhibit density, diffusion and structural 
anomalous behavior as the temperature is increased particles move from the second coor- 
dination shell to the first coordination shell. In the case of systems in which no anomalies 
are present, as the temperature is increased particles will move from the first and second 
shells to further shells. Therefore, in principle the information about the behavior of parti- 
cles up to the second coordination shell would be enough to predict if a system would have 
anomalous behavior. 

(2) 

In order to test this assumption we compute the integrals in the expression for s 2 (see 
Eq. (jSj)) up to the second coordination shell, namely 

s 2 = —^Tcp / [g(r) hig(r) — g{r) + l]r 2 dr, (9) 
Jo 

where r2 is the distance of the second shell. 

Fig. [8] shows the density dependence of S2 along a series of isotherms spanning from 
T* = 0.60 to T* = 1.50 for the cases A 3 , B 3 and C 3 . 

Fig. [9] shows the density dependence of s 2 along a series of isotherms spanning from 
T* = 0.60 to T* = 1.50 for the cases A 3 , B 3 and C 3 . The dots are the simulational data and 
the solid lines are polynomial fits. 

Notice that both s 2 and S2 has a maximum and a minimum for the cases A 3 and B 3 but 
not for the case C3 what indicates an anomalous behavior in the excess entropy. Comparison 
between Fig. [8] and Fig. [9] shows that the excess entropy computed up to the second shell 
not only gives the same trend but also the same density for the maximum and minimum of 
the excess entropy. 

Errington et al. have shown that the density anomaly is given by the condition 
S cx = (<9 s cx/<9hip) T > I 65 - They have also suggested that the diffusion anomaly can 
be predicted by using the empirical Rosenfeld's parametrization? 5 . Based on Rosenfeld's 
scaling parameters and approximating the excess entropy and its derivative by the two body 
contribution, namely, s ex ~ S2 and S ex = S2 the anomalous behavior of the thermodynamic 
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FIG. 8: Excess entropy, s 2 , versus reduced density for fixed temperatures T* = 
0.6, 0.7, 0.8, 0.9, 1.0, 1.1 from the bottom to the top for the case A3, T* = 
0.6, 0.7, 0.8, 0.9, 1.0, 1.1,1.2 from the bottom to the top for the case B% and temperatures 
T* = 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, 1.2 from the bottom to the top for the case C 3 . The 
temperature T* = 0.5 contains densities that are metastable points regarding the high density 
liquid phase. The dot-dashed lines locate the density of maxima s 2 and the dashed, the minima. 

and dynamic quantities are observed if 



excessentropy — > > 

diffusivity ->■ ^ 2 > 0.42 
viscosity — > > 0.83 

density -> ^ > 1.00 . 



(10) 



This sequence of anomalies is consistent with the studies of Yan et al.— iiH&kM. -M where 
structural anomalies are found to precede diffusivity anomalies, which in turn precede density 
anomalies. 

In order to check these criteria in our family of potentials we compute of given by: 

ds 2 



= s 2 -2tt p 2 I In g(r)— — r 2 dr . 
a In p J T J op 



(11) 



And also to test if the criteria given by Eq. (TT0|) can also be applied for computations of 
the excess entropy derivative computed up to the second shell we also calculate 



E 



(2) 



ds 2 
dlnp 



.(2) 



2irp z I ln#(r, 

'0 dp 



d9{r) r 2 dr 



(12) 
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FIG. 9: Excess entropy computed up to the second coordination shell, s|, versus reduced density 
for fixed temperatures T = 0.6, 0.7, 0.8, 0.9, 1.0, 1.1 from the bottom to the top for the case 
A3, T = 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, 1,2 from the bottom to the top for the B3 and temperatures 
T = 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, 1.2 from the bottom to the top for the case C3. The temperature 
T* = 0.5 contains densities that are metastable points regarding the high density liquid phase. The 
dot-dashed lines locate the density of maxima S2 and the dashed, the minima. 

Fig. [TD] and Fig. [TT] show the behavior of and of with the density respectively for 
a fixed temperature for the potentials A3, B3 and C3. The horizontal lines at = 0, 0.42 
and 1.00 indicate the threshold beyond which there are structural, diffusion and density 
anomalies respectively. The graphs confirm that the density, the diffusion and the structural 
anomalous behavior is observed for the potentials A3 and B3 but not for the potential C3, 
confirming Errington's criteria. 

The comparison between Fig. [10] and Fig. [11] shows that the derivative of the excess 

entropy computed up to the second shell is a good approximation for for all the cases. 

(2) 

This result together with the good agreement between S2 and s 2 supports our surmise 
that focusing in the first and second shell behavior we can understand the mechanism that 
leads to the anomalous behavior. 

V. CONCLUSIONS 

In this paper we analyzed three families of potentials characterized by two length scales: 
a shoulder scale and an attractive scale. We found that when approaching the two scales and 
keeping the slope between them fixed, the liquid-liquid critical point goes to lower pressures 
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FIG. 10: YI2 versus reduced density for temperatures 
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FIG. 11: YI2 versus reduced density computed up to the second coordination shell 

while keeping the critical temperature fixed. This result seems to indicate that the slope of 
the curve might be related to the critical temperature while the distance between the scales 
control the critical pressure. This assumption is also supported by another continuous 
spherical symmetric potential in which the slope was varied^. 

We also found anomalous behavior in the density, in the diffusion coefficient, in the 
structural order parameter and in the excess entropy for all the cases in which the distance 
between the scales were not too short. 

From the behavior of the radial distribution function we did infer that the anomalies are 
related to particles moving from one scale to the other. 

In order to check our assumption the excess entropy and its derivative were computed 
in two ways: the total value and the value computed by integrating up to the second 
coordination shell. We found that the behavior obtained by computing these quantities up 
to the second coordination shell is accurate both for system with and without anomalies. 
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Since the Rosenfeld's parametrization 7 ^ is just a lower bound, we can say that it is enough 
to compute S2 and up to the second shell to know if the anomalies are present or not. 
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